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Abstract. We study the magnetoelectric coupling at the surface of a topological 
insulator. We are in particular interested in the surface current induced by a 
static Zeeman/exchange field. This surface current can be related to the orbital 
magnetization of the system. For an insulator with zero Chern number, the orbital 
magnetization is independent of the details at the boundary. With the appearance of 
surface states in the topological insulator, it is not immediately obvious if the response 
is not affected by the conditions at the surface. We investigate this question using exact 
diagonalization to a lattice model. By applying a time-reversal symmetry-breaking 
term near the boundary, no matter if the surface states are gapped out, we still find no 
change in the surface current. This arises from cancelations between Pauli and Van- 
Vleck contributions between surface and bulk scattering states. We also show that the 
surface current response is independent of the chemical potential when it is within the 
bulk gap. Our results are consistent with the claim that orbital magnetization is a 
bulk property. 
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1. Introduction 

The discovery of topological insulator PEIE] has generated strong activities in the 
condensed matter and high-energy physics communities. A prominent property of a 
topological insulator is the existence of gapless surface states in the form of a Dirac 
Hamiltonian with momentum and spin directions locked with each other [H [5]. There 
have been many interesting theoretical predictions, e.g. [61 m El m HO], based on the 
current-spin coupling of this surface Hamiltonian. If the Fermi energy is within the 
bulk gap, it is easy to expect that some physical quantities or phenomena, such as 
the Landau-level spectrum in a perpendicular field [H [HI [I2], transport phenomena 
involving the surface states [13], or the RKKY interactions [HHS] between spins on 
the surface, can be evaluated from an effective surface Hamiltonian. On the other hand, 
one may question if some physical quantities do not just rely on effective surface models 
[16] . The response induced by an external magnetic held is one of the interesting topics 
with above ambiguity. Generally, this perturbation to a system has an effective Zeeman 
held term, a vector potential term, and further contributions [nnHmn]. In this paper, 
we focus on the surface current response to a static Zeeman or exchange held. We 
examine if the surface Hamiltonian is sufficient to determine this response. We also can 
relate the current j to the orbital magnetization MbyM = i/ rx j(r)(ir [20], where r 
is the position. 

Orbital magnetization is an interesting quantity. Even when the maximally 
localized Wannier functions are used, this magnetization, besides a term (“local 
circulation” in [HI [22] or “self-rotation” in [23]) which can be interpreted in the 
same way as rotational motion of electrons in isolated atoms, contains a contribution 
(“itinerant circulation” in [H ] l22 ] l23]) which describes the electronic motion at the edge 
of the sample. Despite the existence of such a current near the sample boundary, it 
has been shown that the total orbital magnetization is independent of the details at 
the boundary for an insulator with zero Chern number [211 E2]. With the presence 
of topological surface states for Chern insulators, this independence becomes even 
less obvious. Nevertheless, Bianco and Resta [2l] provided a real-space expression 
of orbital magnetization for any two-dimensional insulator with finite Chern number, 
and demonstrated that the formula is independent of boundary conditions for a large 
sample. On the other hand, Chen and Lee [25] argued that, for a system composed 
of two insulators with opposite Chern numbers, the orbital magnetization should be 
independent of the details at the boundaries. They also demonstrated, by numerical 
calculations, that the orbital magnetization is insensitive to local perturbed potentials 
near the edge in a Chern insulator, and it is also unchanged with magnetic moments 
at the edges, which gap out the edge states, in spin Hall insulator. Therefore, it is 
interesting to verify if this bulk property still applies for a three-dimensional topological 
insulator, where there are topologically required surface states with the low energy 
physics described by two-dimensional Dirac Fermions. 

In an earlier work, one of us [26] analyzed the surface current of a topological 
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insulator generated by an in-plane static Zeeman field via perturbation theory, 
employing an effective Hamiltonian that is valid only for small bulk momentum k. 
He showed that there is a cancelation between the contribution from redistribution 
of particles (Pauli) and that from virtual transitions (Van-Vleck) between the surface 
states. In particular, due to this cancelation, the current response is independent of the 
chemical potential when it falls within the bulk gap. Furthermore, it is found that the 
bulk scattering states also provides a significant contribution of opposite sign to that 
of the surface states. As a result of this cancelation, the total surface current has no 
dramatic change when the system changes from a trivial to a topological insulator. 

In this paper, to account more carefully contributions from states of all momenta, 
we calculate the current response to a static external Zeeman field, employing a lattice 
model. Using exact diagonalization, we examine the contributions from both the surface 
and the scattering states. To gain more physical insights, we shall consider in some detail 
the contributions from states with small momenta parallel to the surface and compare 
these results with those obtained analytically in |26] . The cancelation between the Pauli 
and the Van-Vleck contributions for the surface states, and the Van-Vleck contributions 
from the surface versus the bulk states, will be examined more carefully. We compare 
the current responses between topological insulators and normal insulators. For the 
topological insulator, we demonstrate the independence of the response on the chemical 
potential (so long as it is within the bulk gap), even though the occupation of the 
surface states depends on it. We shall also show that, due to compensation between the 
surface and scattering states contributions, the total current response is independent 
of surface magnetic moments which gaps out the surface states. We shall also discuss 
effects of different types of external fields on different surfaces. We consider anisotropic 
topological insulators and examine the origin of the anisotropy in the current response. 
Some technical difficulties we encounter in using the exact diagonalization to evaluate 
the Pauli contributions and the contributions near the Dirac point are mentioned. 

This paper is outlined as follows. In Sec. [21 we describe the theoretical model being 
used. The numerical results for the current in response to an in-plane external field are 
provided in Sec. [3l We compare the results obtained here with the analytic study of 
[26] . We also study the effects of chemical potential, anisotropy, and surface magnetic 
moments. In Sec. IH cases with different types of external fields acting on different 
surfaces are discussed. Sec. [5] is the conclusion. 

2. Model Hamiltonian 

An effective model for a topological insulator, for example Bi 2 Se 3 , at small momentum 
k is of the form [1] 

H = moa: + v^k^ay - v{k^Sy - kysfija^. (1) 

Here Si and cXi are the Pauli matrices acting on spin and orbital subspaces, respectively, 
m, V and Vz are material parameters, v and Vz are set to be positive. The parity 
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operator is taken as To include the contribution from finite momenta, we extend 
the Hamiltonian Eq. ([1]) to the following form: 

H = m'{k)ax + sinikfiay - v{sm{kx)sy - sm{ky)sa;)az, (2) 

where m'(fc) = m + c(2 — cos^k^) — cos{ky)) + c^(l — cos(fc^)). For convenience, we set 
the lattice constant (a) equal to one. c and Cz are also material parameters and limited 
to be positive, m > 0 specihes a trivial phase. The condition to be in the topological 
insulator depends on both m, c and c^. We shall focus on the topological non-trivial 
phase within the following region, 

— 2 min{c, c^} < m < 0. (3) 

Here for simplicity, we have written down a form that is appropriate for a crystal of 
tetragonal ^ 4/1 symmetry instead of the symmetry for the real Bi 2 Se 3 family. We 
make this simplication because the relevant issues of the current response in this paper 
is whether the contribution is from the surface states or from the scattering states, 
hence our main conclusions are independent of the lattice symmetry and the specihc 
parameters of the materials. 

We shall start from lattice model to take account of effects of boundaries. The 
corresponding real-space Hamiltonian for a simple cubic lattice is 

H{r) = {m+ 2c +Cz)Y.rCifi^xar 

T Sr [®r ( 2*^® 

Tn^ ( ^ 2 ^x^z)^r+&y 

+af{-^az; - i^(Jy)ar+sz + h.c.], 

where af and are Fermion creation and annihilation operator, respectively, at site 
r. To calculate the surface current, for example, on the top surface with the normal 
along +z direction, we shall apply an open boundary condition to the z direction and 
periodic boundary conditions to the x and y directions. In this case, we get an effective 
Hamiltonian on an one-dimensional chain along the z direction with two good quantum 
numbers, the in-plane momentum components kx and ky. We get the energy spectra 
from the exact diagonalization. For the topological trivial phase, the spectra consist 
simply of scattering states, related to the bulk plane-wave states found in a system 
without boundaries and occupying the same energy range. For the topological non¬ 
trivial phase, there are also mid-gap surface states. These midgap states form a Dirac 
cone for the top {+z) surface and another for the bottom {—z) surface. The energy 
dispersion of the surface states is given by Eg = ±v^Jsirfi{kx) -|- siifi (ky), which can be 
related to zero energy states of its corresponding supersymmetric Hamiltonian. [2^ The 
condition to have surface states is |m|| < jc^l, where = m'{kz = vr/2). Combining 
the condition |m^| < |c| for the side surface in +x direction, where = m'{kx = vr/ 2 ), 
we shall set m within Eq. ([3]) for simplihcation. Note that low energy dispersion at 
half-hlling is described by massless Dirac Fermions. For a given in-plane momentum, 
the decay length of the surface state is A = | In |m|/c^| in the z direction, which 
reduced to n^/|m| as and \m\ Vz- In this paper, we shall focus on the case 
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with the thickness much larger than A so that the coupling between the surface states 
on different surfaces can be ignored. 


3. Current response to an in-plane Zeeman field 


When applying an in-plane Zeeman held By to a system with spin-orbit coupling, one 
can expect redistribution of particles around the Fermi surface because these states 
have spin-momentum locking and the energy levels become lower or higher as the spin 
prefers to be along or against the external held. Therefore there is a contribution to 
the current in the transversal direction of the held. This is the Pauli contribution of 
the current response. However, the Zeeman held also modihes the wavefunctions. The 
modihcation will trigger oh virtual transitions between occupied and empty states |28j . 
This is the Van-Vleck contribution of the current response. As mentioned in Ref. |26] . 
the virtual transition between the surface states is not enough to get a physical answer 
for the current response to an external held, since the result depends on the momentum 
cut-oh employed. A proper evaluation of the Van-Vleck contribution must consider 
virtual transitions from any occupied state to any empty state. Here we use exact 
diagonalization to evaluate all the states and the current they carry. The operator for 
current density along the x direction is given by 


Jx{z, fc||) = dH{z, k{)/dk^, (5) 

where H{z,k\\) is the Hamiltonian after partial Fourier transformation as described 
below Eq. (|1]). Therefore, for the model Eq. (|2]), this is given by csm.{kx)crx — 
V cos{kx)syaz- We obtain the eigenstates with hnite external held By by exact 
diagonalization as mentioned before, and evaluate the expectation value of this current 
operator. For the external perturbation to create the surface current, we use Hsi^y = 
—giSyBy where gi is also a material parameter and would be absorbed in the dehnition 
for By for simplicity. The most general perturbation by an external Zeeman held By 
has two types , with and without an additional ax coefficient in the above expression. 
We shall see that the term with ax contributes no current response to By on the top 
surface. However, it becomes necessary for other cases. We have more discussion about 
this point in Sec. 01 We shall consider systems with a uniform By at hrst, and discuss 
the ehects of non-uniform By in Sec. 13.41 For a uniform By, the current deep in the 
bulk vanishes, since Eq. (0]) for the bulk has inversion symmetry. The current hows 
only near the surface. The surface (number) current on the top surface can then be 
evaluated by 


F = 

^top 


r d% 

Ib.z. 


jGtop half 


( 6 ) 


where k\\ is the momentum parallel to the surface. The surface current for the bottom 
surface is from similar formula but the sum z is over the bottom half. Although we 
could have studied general responses, we shall conhne ourselves to the linear case where 
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Figure 1. (Color online) Current response of a topological insulator from states with 
fc|| < kr- (a) Total contributions, (b) Contributions from the bulk scattering states, 
(c) Contributions from the surface states, c = v = Cz = Vz = ^- Kb’s are the same 
for different g’s. A sharp change of Kg around the Fermi energy is corresponding to 
the Pauli contribution. We limit the chemical potential to be 0 < /i < |m|, where 
m = —0.4 in these plots. The magenta dash line is linear fitting around small k for 
the case with g = 0. 


physical pictures can be accessed more easily. We therefore present our results in the 
form of a linear response coefficient dehned by 

K = UJBr m 

3.1. Surface and scattering states contributions 

3.1.1. Chemical potential: First of all, we shall show that the current response to HBi,y 
is independent of the chemical potential, provided that it is within the bulk gap. In 
this subsection, the chemical potential independence is verihed from numerical results. 
Combining with other studies below, we shall also give a physical argument why this 
response is independent of the chemical potential before the end of Sec. 13.31 

In order to make the physics more transparent, we shall hrst compare the results 
from numerical calculations with those from analytic works [26] by studying the response 
of the states with small k\\. We do this by limiting the domain in Eq. (j6|) to k\\ < kr'. 

[ d%i 271 dh 


B.Z. 


0 


( 8 ) 
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where k\\ = |A;||| is the magnitude of the in-plane momenta. We divide the total current 
response {k) into two contributions, one from the surface states (denoted by Kg) and 
one from the bulk bands (denoted by Kb). The results are shown in Fig. [H We limit the 
chemical potential to be |/i| < \m\, so that the system has no scattering states around 
the Fermi energy. Figure [11(b) shows that Kb is independent of /i, though it is hnite. 
This independence is in accordance with [26], where it was concluded that there is no 
contribution to the current via virtual transitions from the valence band to the surface 
states and from the surface states to the conduction band. That is, Kb arises entirely 
from virtual transitions between the valence and conduction bands. Kg are shown in 
Fig. HKc). For /i = 0, there is no Pauli contribution due to the vanishing of density 
of states, and Kg is entirely due to virtual transitions from occupied to empty surface 
states. The magnitude of Kg increases with kr when more surface states are included, 
as shown by the dashed line (Our actual numerical calculation gives the yellow dots 
with some small deviations from this line near small kj.. This is related to a technical 
issue which will be discussed in Sec. Id.l.dh . For hnite p’s. Kg becomes zero when kr 
is less than the Fermi momentum kp. This rehects the absence of virtual transitions 
between the surface states when kr < kp due to the Pauli exclusion principle, since at 
a given momentum k\\ with k\\ < kp, both the surface states with Eg^O are occupied. 
Around kr ~ kp, the response has a sharp change, which is mainly from redistribution 
of particles by the external held. This represents the Pauli contribution. The hnite 
width of the change in Kg in Fig. [1] is due both to the fact that the Fermi surface is 
not a regular circle and the hnite the resolution of the points that we have used. This 
Pauli contribution complements exactly the missing Van-Vleck contributions (compared 
to p = 0) from surface states with k\\ < kp. Kg becomes independent of p when kr > kp. 
K = Kg + Kb given in Fig. [I](a). At these small momenta. Kg is main contribution to 
K, while Kb is a smaller contribution of opposite sign. However, as shown in Fig. [3l Kb 
shows non-monotonic behavior and has comparable amplitude to that of Kg for some 
cases. The comparison between Fig. [1] and Fig. [3] for isotropic case also shows that 
Kb versus Kg diher more (less) for increasing (decreasing) |m|. The numerical results 
for small kr discussed above is consistent with the analytic results in [26]. The current 
response is insensitive to the occupancy of the surface states, which infers it is a quantity 
related to the bulk properties. 

Numerical results also show that k for the bottom surface has equal amount but 
opposite sign to the top surface. This is related to the charge conservation. We shall 
provide an argument in Sec. 13.31 why this current must be independent of the chemical 
potential. [|| 

f An alternate way to understand this result is as follows. As mentioned in Sec. 13.21 we can consider 
our system as a collection of two-dimensional systems, one for each ky (for a magnetic field along y). 
In our case, each of these two-dimensional systems has zero Chern number (at finite By, as can most 
easily be seen by considering the surface state spectra). On the other hand, the derivative of the orbital 
magnetic moment with respect to the chemical potential is proportional to the Chern number for a 
general two-dimensional system |22j . Hence this derivative must also vanish for our three-dimensional 
system. 
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Figure 2. (Color online) Current response for different m’s to By on the top surface, 
(a) Fixed c and varying c^. (b) Fixed Cz and varying c. Total current responses on the 
side surface are found to be the same as those on the top surface. 


3.1.2. Dependence on Hamiltonian parameters: The total surface current responses 
are dependent on the material parameters, i.e. m, c, c^, v and Vz- We like to show 
more about this point. For convenience, we set c = v and = Vz- Because of the 
independence of k on the chemical potential, we shall study the systems with /i = 0 . 
For a system with c 7 ^ Cz, it means the anisotropy between in-plane and out-of-plane. In 
Fig. [2], we show that the system has different current responses and anisotropic effects 
between topological and normal insulators. In the normal phase, non-zero current from 
only the scattering states is slightly modihed as changing the material parameters. As c 
or Cz become larger, the responses are larger. When m changes from positive to negative 
values, the response continually grows up into the topological phase. Instead of weak 
dependence and monotonic increasing for the current response of normal insulator, it 
shows different dependence on c and in the topological phase. Figure [3] shows more 
details about the anisotropic effects. Note that the external held is set to be so small 
that the surface states are well separated from the bulk bands for all momenta. In the 
normal phase, the current response is only from the scattering states and enhanced for 
larger c or c^, as shown in the middle row of Fig. [3l The amplitude of the response is 
monotonic increasing as a function of kr for small momenta. For the topological phase, 
the dependence is more complicated. First, we focus on Kg in the bottom row of Fig. 
[3l Before saturation, the slopes of Kg are the same for all cases, which related to the 
independence of c for the Van-Vleck contribution. m We also hnd that total Kg is 
weakly changed as modifying c^, but strongly dependent on c. This might be related to 
the fact that the size of surface cone becomes smaller as c is larger. However, the results 
of K (top row of Fig. [3]) show that it is still dependent on c^, which mainly related to kj, 
(middle row of Fig. [3]). Kf, shows an opposite contribution for those momenta with the 
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Figure 3. (Color online) Anisotropy study for topological (black line) and trivial 
insulator (broken red line) on the top surface. The external field is along the y direction. 
Three cases with different material parameters are shown in different columns. Top row 
is for total current response, middle row for the contribution from scattering states and 
the bottom for that from surface states. Black full lines are for topological non-trivial 
phase and red dash lines for trivial phase. 


surface states. Once the surface states merge into bulk bands, Kb is changing in a way 
similar to that of normal insulator. We note that k increases smoothly as crossing the 
transition, which infers that Kb compensates the missing of the midgap states at larger 
momenta. Basically, the momentum region to have surface cones is a bulk property. 
In addition, k is independent of the chemical potential as |/x| < \m\. Therefore, the 
current response to an external held cannot be described only from effective surface 
Hamiltonian. 

3.1.3. Technical remarks: Before further discussions, we like to mention some technical 
problems we encounter when calculating the current by the exact diagonalization, and 
how they are overcome. The hrst is the difficulty to calculate the contribution with 
momenta near the Dirac point as /i = 0. Due to hnite size effect, two surfaces couple 
to each other and result in a gap, inversely exponentially proportional to the size of 
the sample, to the ideal Dirac cones. This causes a sudden jump for the current as the 
energy scale of the external held is larger than this gap. The corresponding current 
response will rehect size ehects and this is already beyond the linear response, which 
is not the concern in this paper. From a physical point of view, the density of state 
approaches zero at the Dirac point, and so the Pauli contribution should vanish. In 
the numerical calculations, we either subtract out the spurious contributions to the 
current (when larger helds were used in the exact diagonalization), or slightly shift the 
momentum to avoid the ambiguity. Another problem is the difficulty to count the Pauli 
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Figure 4. (Color online) fcy-resolved current response for the top and side surface. 
The current has the same value for different surfaces for each ky. c = 1, Cj = 1.3, 
m = —0.4. 


contribution for cases with p ^ 0. Finite grid makes the evaluation of redistribution of 
particles difficult. We can overcome this problem by applying twist boundary conditions 
and averaging over different cases. This method is equivalent to including more lattice 
momentum points. The twist boundary condition can also be applied to estimate the 
contribution around the Dirac point. 

3.2. Relation to 2D system 

Now, we study the current response on different surfaces. We shall focus on the —x 
surface as a side surface {i. e. surface normal pointing from the sample to its outside along 
—X ). The current generated by By on this side surface is along the — 2 ; direction 
and calculated from a formula analogous to Eq. (|6]) with x and z interchanged. In this 
case, we apply an open boundary condition along the x direction and periodic boundary 
conditions in the y and z direction. Though the calculations in this subsection is done 
for a him geometry which is different from those in the last (Sec. 13.ip . our calculations 
show that the current response is identical for the same given parameters c and (Fig. 
Ep. This is expected, since for a sample with hnite size in both x and z directions, 
we must have a (number) current circulating counterclockwise along the y direction. 
The currents on different surfaces must be equal, which is inevitable due to the charge 
conservation. Even more we hnd, as shown in Fig. 01 the current is conserved for each ky. 
{ri{ky) is dehned as by Eq. ([7]) but with the numerator replaced by Ifopi^y) or Igideiky), 
which are in turn given by a formula similar to Eq. (E]) except there is no integration 
over ky.) It is easy to get this conclusion for ky = 0, because the system reduces to a 
two-dimensional quantum spin-hall system. For hnite ky, one has a tight-binding model 
for a hctitious two-dimensional system with time-reversal-symmetry breaking terms. 
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However, the currents must still be the same for different surfaces again due to the 
charge conservation. Our results give a numerical demonstration that current for each 
momentum ky is conserved. This conclusion applies to the cases with surface magnetic 
moments, which is discussed further in the next subsection. 

Figure m also shows that the current response is mainly from small ky and decays 
as increasing \ky\. We note that K{ky) is a smooth function in ky, even though surface 
states exist only for \ky\ below some maximum value (where the surface states merge 
into the bulk states). The sum of the surface and scattering states contributions give 
rise to smooth ky dependence. 

The current responses of the side surface for the cases with different m and 
anisotropy are also the same as those of the top surface for each given ky. The resulting 
k’s are the same between different surfaces, as mentioned in the caption of Fig. |2l 
However, the dependence of the current response on kr is different. Numerical results 
for small momenta show that k on the side surface is weakly enhanced by increasing 
and strongly suppressed by increasing c. This might be explained by the shrink of the 
momentum region to have the surface states and the change of slopes of the surface cones. 
However, this momentum region is related to bulk properties of a topological insulator 
by material parameters, as discussed in the last subsection. In addition, anisotropic 
effects of the current response to on side surfaces are not simply explained by the 
modihcation of the surface cones. We have more discussion about this point in Sec. ID 

3.3. Surface magnetic moment 

Even with the independence of the chemical potential, one may question if the current 
changes dramatically when the surface states are gapped out, for example in the presence 
of an extra term in the Hamiltonian proportional to s^, representing local moments, at 
the top surface. For clearer explanation, we shall focus on the isotropic case (c = = 1) 

and study the current response to By on the top surface. The responses to other helds 
on different surfaces can be found by a mapping discussed in the next section. 

Before calculation, it is helpful to apply some simple arguments to hnd out what 
we should expect. Consider surface magnetic moments put just on the top layer. If the 
size of the system is large enough in the 2 ; direction, the surface current for the bottom 
would not be affected. Our numerical study has verihed this point for all considered 
cases of this paper. Using charge conservation, the current at the top surface cannot 
change either. (This point may be clearer if one considers a sample with hnite (but 
both large) sizes in both the x and z directions, since, as mentioned in the last section, 
we must have a circulating current). Hence, the current at the top surface cannot be 
affected by the surface magnetic moment, even though a gap may appear. (The same 
conclusion can be reached by consider the currents at the side surfaces.) We shall verify 
this below by explicit calculations. 

To include local magnetic moments at the surface, we add an additional term 
Hm = hi(r) ■ s 


( 9 ) 
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Figure 5. (Color online) Current response for the system with surface magnetic 
moment hi. (a) The robustness of the current to local magnetic moments for each ky. 
The contributions from the scattering and the surface states are shown in (b) and (c). 
Same symbols are used for all panels. The blue dash line in (c) is linear fitting around 
small k for the case without surface magnetic moment and with g = 0. m = —0.4. 


to the Hamiltonian, where hi(r) is an position dependent exchange held. (This is not 
the most general form (see Sec. 0] below) but is sufficient for our considerations here). 
We hrst consider an hi(r) to produce a gap on the surface states. For this we use 
hi(r) = hi{r)z with hi uniform and hnite only on the top layer of our system. (The gap 
generated is however smaller than hi by a factor from the overlap between the surface 
state wavefunction with the hrst layer). Numerical results of the current responses for 
various values of hi are shown in Fig. |5l Figure [5](a) shows that the total current 
response is unchanged: in fact, the /Cj^-resolved current is also independent of hi. (Note 
that while Fig. 0] presents the current conservation between diherent surfaces. Fig. [5](a) 
shows the robustness of the current to surface magnetic moments.) 

The above however does not imply that the current contribution from each in-plane 
momentum unahected by hi. We show the kr dependence at small momenta in Fig. 
|5](b-c). First, let us focus on /i = 0. With hnite hi, Kg is suppressed. This can be 
easily understood by considering the modihcations to the surface state spectrum when 
a gap is opened up. On the other hand, Kf, also decreases, and in fact can even change 
sign at small kr. These changes compensate each other exactly when integrated over 
all momenta, giving a response independent of hi. Figure [5] also contains a comparison 
between /i = 0 and fi = 0.2 at a hnite hi. As in Sec. 13.11 Kg vanishes exactly for 
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kr below the Fermi momentum, implying once more that the Van Vleck contributions 
involving virtual transitions between the surface and bulk bands vanish even in the 
present case, even though the boundary condition at the surface becomes different from 
that of Sec. 13.11 since now there is a hnite hi at the surface. When is increased 
beyond the Fermi momentum, the Pauli contribution again compensates exactly the 
missing surface-to-surface Van Vleck contribution, as in Sec. 13.11 

From the physical argument given in the beginning of this subsection (or the direct 
evaluations above), the independence of the current on the chemical potential in Sec. 
13.11 is then no surprise. Consider again a given sample with magnetic moments on the 
top surface, gapping out the surface states. The current at the top surface should be 
independent of /i if /i is within the created gap. By charge conservation, the same must 
then be true also for the ’’pristine” bottom surface. One can consider the large hi limit 
where the surface states disappear entirely for the top surface and conclude that the 
total current at the bottom (and hence for all surfaces) must also be independent of fi. 
The same arguments can be applied to the side surfaces when considering a sample with 
hnite (but large) extent in both x and z directions. In this case we must have circulating 
current of equal magnitudes near both the x and z surfaces when the magnetic held is 
along the y direction at all /r’s within the gap and independent of whether the surface 
states on one particular surface is destroyed or not. 

3 . 4 . Non-uniform external fields 

This robustness of the current response to a time-reversal-broken term also applies to 
other situations. We shall discuss two special cases, both related to distributed external 
helds. The hrst is an artihcial case with an external held on all sites except the top 
layer. The results are shown in Fig. |6](a-c) by lines with crosses. This case can also be 
considered a system with surface magnetic moment with hi hnite only on the top layer 
and along the —y direction with an magnitude equal to that of the external Zeeman held 
By. Therefore, it is instructive to compare this case with Sec. 13.31 Comparing Fig. [6](b- 
c) with Fig. [5](b-c), we indeed hnd similar dependence of the responses on kr- At small 
kr, Kg is suppressed, but Hb also decreases, which in fact can also change sign near small 
kr- Their sum k at small kr is smaller, but the loss at small momenta is complemented 
by the contribution from large momenta, which mainly from the scattering states. After 
summing over all states, the total current is not ahected as shown in Fig. EKa). 

The second case is a system with By just on the top surface. Because one layer is 
very small as comparing to the size of the bulk, we can expect no total surface current 
(since the current at the bottom surface is necessarily zero). Indeed as shown in Fig. 
Ela), we obtain vanished k after summing over all states. However, when fey is summed 
only up to a hnite range kr, it is in general hnite. Kg at hnite kr is similar to the 
uniform By case, being always negative but smaller. It is somewhat difficult to separate 

§ The same argument applies to the current response to an external orbital magnetic perturbation. 
Therefore we believe the total current response (orbital magnetization) is a bulk property. 
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Figure 6. (Color online) Current response for the system with position dependent 
By. (a) Total current response for the system with external field for all sites (black full 
circle), for all except the top layers (green cross), and just for the top layers (red empty 
diamond). The corresponding contributions from the scattering and the surface states 
for momenta up to kr are shown in (b) and (c), respectively, (d) Position dependence 
of current response with different mass to By just on the top layer, m = —0.4 for 
(a)-(c). 


clearly the surface versus bulk states when ky becomes comparable to the k\\ where these 
states merge into each other, but our calculation shows that Kg seems to vanish when 
summing over all surface states. We understand this as follows. We note that at the 
parallel momentum where a surface merges with the bulk states, the decaying length 
of the surface state diverges, and therefore a finite localized magnetic moment at the 
surface will result only in a vanishing perturbation. The vanishing of the current from 
the surface states can then be understood as due to the fact that now the surface state 
spectra simply shifts in momentum space with their end points hxed. (Note that this 
shift in momentum is therefore k dependent and is different from that from a uniform 
gauge field). From Fig. [6l within ky is also in general finite. In particular, it acquires 
a negative value near the momentum where the surface states merge into the bulk, 
though the total k integrated over all momenta again vanishes. 

Since we are evaluating the linear response, we can see from Fig. [6]that the response 
for a uniform By is just the sum of the responses for the two non-uniform By cases 
discussed above, as it should be. We have also checked that, for a finite p, there is no 
response for ky smaller than the Fermi momentum (as in Fig. [T](c) and[5](c)), suggesting 
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Table 1. Current response of a topological insulator and reaction of its surface states 
to different types of external field on different surfaces. Sh±i: surface cone shift in 
±/ci direction. SI: Slope reduced for all momentum. Sle^: Slope reduced except along 
kx = 0. Sley: Slope reduced except along ky = 0. 


component 

current response 

reaction of 

component 

reaction of 

on +z surface 

on +z surface 

top surface states 

on +x surface 

side surface states 


yes 

Sh y 

Sy 

Sh, 

Sy 

yes 

Sh, 

^Z^X 

Shy 

Sz 

no 

Gap 

^X^X 

Gap 


no 

Sley 

^y^ X 

Sle. 


no 

Sle. 

-Sz 

Sley 


no 

SI 

Sx 

SI 


again that there are no Van-Vleck contribntions dne to snrface to bnlk virtnal transitions 
for onr non-nniform By^s. 

Fignre [6]^d) shows the position dependence of the total cnrrent response for the 
case where By is localized at the snrface. There is a sharp (lattice-scale) featnre near 
the snrface. This is followed by a slower feature of opposite sign decaying into the bulk, 
with an m-dependent decaying length (around 2a for m = —0.4 and 6a for m = —0.15.) 
This pattern then is roughly a current loop around the localized external held By. For 
a thick sample, there is thus no net current on the bottom surface. 

4. Types of perturbed Hamiltonian for external fields on different surfaces 

When comparing the current responses to external helds with different directions, we 
should take more care about the form of the perturbation Hamiltonian. Actually, Eq. 
([2]) or Eq. (01) already implied the anisotropy of topological insulator between the z 
and in-plane directions. For this system, indeed the current response to By is from 
HBi,y If the applied held is along the 2 ; direction, however, the current response of a 
side surface is not simply from SzB^ with some coefficient. In this section, we shall have 
more discussion on this point. 

We can limit possible terms in the Hamiltonian by symmetry arguments. Since an 
external magnetic held is a pseudo-vector which also breaks the time-reversal symmetry, 
the allowed terms in the Hamiltonian describing the ehects of a magnetic held can be 

a 

Hb = —bi ■ s — b 2 ■ scTx, (10) 

where bi and b 2 are external magnetic helds scaled with appropriate directional 
dependent ^f-factors. We shall call the term in Eq. flTOj) without parity operator Hbi 
and the other Hb 2- 
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Our numerical calculations show that, for the top surface, current can be generated 
from bx and by by s^ and Sy, respectively. On the other hand, non-zero current response 
on the x-side surfaces is from by and bz by Sy and SzCx, respectively. These conclusions 
agreed with those of [26] which analyzes only states with small momenta. The responses 
might be possibly be understood by considering the reaction of the surface states to the 
external helds, which is summarized in Table [H Therefore, it might be reasonable to 
expect that the system has a current response as the surface cone is shifted by external 
helds. For example, the component Sy shifts the surface cone in +kx direction and creates 
a contribution of the number current in —x direction if the perturbed Hamiltonian is 
Hsi^y Also, Sy shifts the surface cone in —kz direction for +x surface. Therefore 
Hsi^y will create a surface current along -I- 2 ; direction on the x-side surface, as discussed 
before. Those denoted by gap or reduced slope have no current response. For example, 
another possible term from By is HB 2 ,y, which excite no current on both the top and 
side surfaces. Even though for the total physical response one must consider both the 
surface and bulk states, it turns out that the above criteria can decide whether there is 
or is not a hnite response to a particular external perturbation. 

As a comparison, we like to demonstrate anisotropic effects of the current responses 
to Bz on side surfaces. The maximum momentum region to have surface states within 
the constraint Eq. (]3|) is similar to an ellipse as c 7 ^ c^. The radius of this region in the kz 
direction is set to be pz and that in the ky direction py. From our numerical calculation, 
we hnd that, for systems with hxed c, pz is shorter while Cz is larger. For systems with 
hxed Cz, increasing c will decrease py. The parameter change helps to understand the 
current response of systems with these material parameters. The hrst and third columns 
of Fig. [7] show that k decreases as Cz becomes larger. This might be explained by the 
shrink of pz, which supposed to reduce Kg- Kg no longer contributes to k when is 
larger than some critical value, but the curve of k smoothly increases as crossing the 
transition. Kb compensates the loss in the region without those midgap states. On the 
other hand, the hrst two columns of Fig. [7] show that k becomes larger when c is larger. 
One might relate this increasing to sharper slope of the surface cone in ky. However, the 
shrink of py would make the surface states have less region to contribute Kg. The same 
reduction occurs for Kb. After summation, k still smoothly increases as including more 
contribution from larger momenta. The current responses on a side surface cannot be 
simply explained from the study of the surface states. Besides, the momentum region to 
have surface states is determined by the material parameters, which are used to describe 
the bulk properties. It is necessary to consider both the surface and the scattering states 
in current responses. 

5. Conclusion 

In conclusion, we studied the current response of a topological insulator to a static 
Zeeman held using an exact diagonalization to a lattice model. The effects of current 
responses to different types of external helds on diherent surfaces are discussed. We 
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Figure 7. (Color online) Anisotropy study of the current responses to on a side 
surface. The arrangement of these plots is the same as that illustrated in Fig. [3]except 
the external field is along the 2 direction. 


find no change in the surface current by changing the occupancy of the surface states or 
applying a time-reversal symmetry-breaking terms near the boundary. These only affect 
the individual contributions to the current, such as Pauli versus Van-Vleck, or surface 
versus bulk states, but not the total. This suggests that there is a kind of sum-rule for 
this response, though we have not yet been able to derive it analytically. [JJ] Our results 
support the conclusion that magnetization is a bulk property, independent of the details 
at the boundaries of the sample, even for topological insulators. 
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